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Melting of an ultrathin lubricant film confined between two atomically flat surfaces is we 
studied using the rheological model for viscoelastic matter approximation. Phase diagram with 
domains, corresponding to sliding, dry, and two types of stick — slip friction regimes has been 



built taking into account additive noises of stress, strain, and temperature of the lubricant, 
-t— > 

The stress time series have been obtained for all regimes of friction using the Stratonovich 

interpretation. It has been shown that self-similar regime of lubricant melting is observed when 

intensity of temperature noise is much larger than intensities of strain and stress noises. This 

regime is defined by homogenous distribution, at which characteristic stress scale is absent. We 

study stress time series obtained for all friction regimes using multifractal detrended fluctuation 

analysis. It has been shown that multifractality of these series is caused by different correlations 
> 

that are present in the system and also by a power-law distribution. Since the power-law 
re- 
distribution is related to small stresses, this case corresponds to self-similar solid-like lubricant. 
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1 Introduction 

Problem of sliding friction is of great interest due to its applied engineering significance 111 . Atom- 
ically flat surfaces separated by ultrathin layer of lubricant are under active investigation recently. 



These systems exhibit anomalous behavior, consisting in existence of several kinetic regimes of fric- 
tion. Transitions between the regimes are interpreted as a phase transitions |2|. The liquid lubricant 
shows properties of solids |3|. Distinctive peculiarity of the systems is inherent in dry friction inter- 
rupted motion (stick — slip regime) (4}j8]. Denoted regime is observed for lubricant thickness equal 
or less than three molecular layers, and is explained by periodical solidification due to walls press- 
ing. Sheared lubricant melts when shear stresses a are larger than the critical value a c (yield point) 
owing to "shear melting" effect. The increased interest to such systems has motivated appearance 
of several models. Deterministic model |6], thermodynamic model |9|, and rheological model (To) 
were developed to describe above mentioned properties. Investigations are also based on molecular 



dynamics methods 11 . The influence of additive non-correlated noises of basic parameters 12, 13 



and correlated fluctuations of temperature fl4) on lubricant melting has been investigated within 



the framework of rheological model flQ| . Reasons for hysteresis behavior 115-17 and melting due 
to dissipative heating of friction surfaces [18] have been also considered. Systems with different 



viscosity dependence on temperature are also analyzed 19 



We suppose that with the increase in stress a the lubricant melts, since the velocity of moving 



surfaces also increases according to the relationship 18| 20 



Veff 



where h is the thickness of lubricant or distance between friction surfaces, i] e ff is the effective 
viscosity, being measured experimentally |20|. 

The present work is devoted to time dependencies investigation of stresses in a self-similar regime 



of lubricant melting, caused by temperature fluctuations. This regime was found in 13 based on 



the method described in 21 , 22 



2 Dynamic phase diagram 



In previous works 10, 12, 13, 18 we treated a viscoelastic medium with a non-zero thermal con- 



ductivity using the rheological model. The system of kinetic equations was also derived describing 
mutually coordinated evolution of shear stress a and strain e, and temperature T in ultrathin 



lubricant film during friction between atomically flat solid surfaces. We used the measure units 

( pc v -q T c \ 1/2 a s 

°a= , £ s=7T> T c ( 2 ) 

\ T T J G 

for variables er, e, T, respectively, where p is the lubricant density, c v is the specific heat capacity, 
T c is the critical temperature, rjo = i] at T = 2T C is the characteristic value of shear viscosity rj, 
Tt = ph 2 c v /n and h are the time of heat conduction and thickness of lubricant, k is the thermal 
conductivity, r £ ~ 10~ 12 s is the relaxation time of matter strain, Go = ?7o/ r e is the characteristic 
value of shear modulus. Let us write the equations: 



T <T &=-o a + ge + y/I ff £ 1 (t), (3) 

r e e = -e + (T-iy°+V^2(<), (4) 

r T f = (T e - T) - a a e + a 2 " + V^ 3 (t). (5) 



Here the stress relaxation time r ff , the temperature T e of atomically flat solid friction surfaces, and 
the constant g = G/Gq < 1 are introduced, where G is the lubricant shear modulus. Quantities I a , 
I £ , and It are the intensities of stress, strain, and temperature noises, respectively. Substitution 
of de/dt instead of e/r a in Eq. p| reduces it to a Maxwell-type equation for a viscoelastic matter, 
which is widely used in the theory of boundary friction 1 1 1. The relaxation behavior of a viscoelastic 
lubricant during friction is also described by the Kelvin- Voigt equation (El) 10 23). It takes into 



account the dependence of the shear viscosity on the dimensionless temperature r] = f]o/{T — 1). 
Equation ^ is heat conduction expression describing heat transfer from friction surfaces to the 
layer of lubricant, the dissipative heating of the stress- induced viscous flow, and a heat source due 
to the reversible mechanocaloric effect. Equations (J3j) - (J5J) formally coincide with the Lorenz 



synergetic system 24,251, where the shear stress acts as the order parameter, the conjugate field is 



reduced to the shear strain, and temperature is the control parameter. When a = the lubricant 



is solid-like, situation with a ^ corresponds to its liquid-like state 10, 13-19]. 

In equations (J3J) - ^ 0<a<l is the fractional exponent. The function £j(£) is ^-correlated 
Gaussian source (white noise). Its moments are defined aaM 

(&(£)> = 0, utMt')) = ^At-t'). (6) 



1 Here multiplier 2 is chosen for simplification of the corresponding Fokker - Planck equation (FPE). 



Experimental data for organic lubricant |3| show that relaxation time of the stress r CT at normal 
pressure is ~ 10~ 10 s, and it increases by several orders of magnitude at large pressures. Since the 
ultrathin lubricant film consists of less than three molecular layers the relaxation process of the 
temperature to the value T e occurs during time satisfying condition r^ <C r a . Then, within the 
adiabatic approximation r CT ^> r e , tt, equations (El) and (J5J) are reduced to the time dependencies 



e(t) 



e + eUt), T(t)=T + T&{t); 



(7) 



a a (T e -l+a 2a ) d a (a), s=Vl £ +I T a 2a d a (a 



f = (T e +2a 2a ) d a (a), T = ^I T +I e a^ d a (a), d a (a) = (1 + a 2a )-\ (8) 

Here, deterministic components are reduced to expressions obtained in Ref. |10|, whereas fluctua- 
tional ones follow from the property known as variance additivity of independent Gaussian random 



quantities 26 . Thus, the use of the slaving principle inherent in synergetics |25] transforms ini- 



tially adiabatic noises both of strain e and temperature T to multiplicative form. As a result, a 
combination of Eqs. (J3j), ([7]), and ([8j leads to the Langevin equation 12 , 13 , 18 : 



& = f a (a) + ^Ua)at), 



(9) 



where the time t is measured in the units of stress relaxation time T a . Generalized force f a (a) and 
effective intensity of noise I a {c) are fixed by equations fl2J|T3~|[l8~] : 



(10) 



f a (a) = -a a + ga a [l - (2 - T e )(l + a 2a y l ] , 
h{°) = I a + g 2 {I £ + I T v 2a ){l + v 2a y 2 - 



Effective intensity of noise is obtained in accordance with variance additivity property of noise 
mentioned above. In order to avoid mistakes, one should notice that a direct insertion of Eqs. ^\ 
and (pi) into (J3J) results in the appearance of a stochastic addition 



l¥ 2 + 4 1/2 + 4 7 V )gd a (a) m 



:n) 



whose squared amplitude is quite different from the effective noise intensity (10). Moreover, in 



contrast to the expressions (J8l), a direct use of the adiabatic approximation in Eqs. Q and (J5J) 
reduces the fluctuational additions in Eqs. (TH) to the forms: e = [I e + I T cr a )d a (a), T = (I T — 
I £ a a )d a (a). The latter is obviously erroneous since the effective noise of the temperature T 



disappears entirely for the stress a = (7r// e ) 1 / 2a . The reason for such a contradiction is caused by 
the fact that Langevin equation does not permit the use of usual analysis methods (see |26|). 

Langevin equation (J9J) is a stochastic differential equation (SDE), since it contains stochastic 
force yJl a (o) £(£). Therefore each solution of the equation is individual and we can say only about 
statistical characteristics of such solutions. In this context further we consider only probability 
distribution of solutions P a (cr) over stress value a. 

Multiplying (J9J) by dt, the Langevin differential relationship is obtained: 



da = f a (a)dt + y/Wa)dW(t), (12) 

where dW(t) = W(t + dt) — W(t) = £(t)dt is the Wiener process with properties J27J: 

(dW(t)} = 0; ((dW(t)) 2 } = 2dt. (13) 



In general case, infinite number of the FPE's forms can correspond to equation (12). 

There are several forms of interpretation, each can be characteristic of specific physical object. 
Three the most used are the Ito interpretation (I-form), the Stratonovich interpretation (S-form) (25) 



and the kinetic form (K-form) |28|. Within the framework of the Ito form stochastic processes a(t) 



and dW(t) presented in the last term of equation (12) are supposed to be statistically independent 



25|. Integrating (12) with the use of the Stratonovich interpretation it is necessary to evaluate the 



last term in the center of the time interval, i.e., to use the following construction 25,27 




dW{U). (14) 



The appearance of correlation between processes a and dW can be seen in this case suggesting 
about the presence of memory effects. Such effects are often present in real systems. In other words 



the Stratonovich form corresponds to the equation (12) with real noise, that can be approximated 
by the Gaussian white noise. In general case the Ito form is used for biological systems with 
discrete time |30] , for statistical interpretation of birth-death processes of living organisms for 
example. Thus, for description of melting of ultrathin lubricant film, one would rather choose the 
Stratonovich form, since in this case the time is continuous and correlations are present in the noise. 



In the works 13 , 18 the simple Ito form of the FPE was used. Here within the framework of the 



Stratonovich form, we show that for our system the change of the form of interpretation does not 



lead to qualitative changes of its behavior. The corresponding FPE with respect to (13) is (S-form): 






d 



y/lJjr)j^y/h(<r)Pa{(r,t) 



(15) 



Distribution of the solutions of (12) becomes stationary in time, and its form can be found from 
g5) at dP a (a,t)/dt = 0: 

P a (cx)=Z- 1 exp{-U a (a)}. (16) 

The obtained distribution is defined by a normalization constant Z and an effective potential 



UJa) 



1 



lnl„ 



<j 



faW) 



da' . 



(17) 



Extremum points of the distribution (16) are defined by condition dU a /da = dl a /da — 2f a = 0, or 
in explicit form 



T e -2 
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So, extremum abscissas of P a (cr) are independent of noise intensity I a . Expression (18) differs 



from analogous one, obtained in [13). In [13] the second term is multiplied by 2a, but in this case 
it is multiplied by a. Thus, at increase in all noises intensities by two times further examination^ 
within the framework of the Stratonovich interpretation concurs with results, obtained using the Ito 
form jl3]. However, potential (17) does not take earlier obtained form [13] at simple renormalization 
of the noise intensities, since it differs from above only by the first term (presence of factor 1/2). 
Therefore the time dependencies of the stresses are different. Since the aim of this work is to study 
peculiarities of time evolution of the stress, we use the Stratonovich approach. Earlier studies [13] 
were focused only on stationary states using the Ito approach. The typical phase diagram at fixed 
temperature T e is shown in fig. [TJ where lines correspond to the stability loss limits of the system. 
Straight line going from the beginning of coordinates is defined by 



2I £ . 



(19) 



This follows from (18) and limits existence of zero stationary solution ctq = 0. Above this line the 



maximum of P a {<?) always exists at o"o = 0, below it the maximum does not exist. In the diagram 
four domains with different regimes of friction can be seen. 

Unnormalized probability distribution ( 16 ) shown in fig. [2] corresponds to domains in fig. fl| Point 



2 We study the extremums of distributions for phase diagrams analysis and interpretation of the stationary states. 
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Figure 1: Phase diagram at g = 0.9, T e = 1.5, a = 0.75 with domains of friction such as sliding 
(SF), dry (DF), stick - slip (SS), metastable and stable sliding (MSF + SF). 
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Figure 2: Distribution ( 16) at I a = 10 10 and regimes shown by points in fig. QJ 1 — I £ — 0, It — 4 
(DF); 2 — I £ = 0,I T = 15 (SS); 3 - I e = 2,I T = 2 (£F). The insert shows P a (a) at 7 e = 
1.968, I T = 3.5 (MSF + SF). 



1 is located in dry friction region (DF) of the phase diagram, and single maximum of distribution 
is observed at o"o = 0. Two-phase region SS is defined by existence of distribution maximums 
P a (a) at zero and non-zero stress values (point 2). Point 3 is located in domain, where only one 
distribution maximum exists at o~ 7^ 0, being related to liquid friction or sliding regime (SF). In 
the insert distribution in area MSF + SF is shown. Dependence P a (a) has two maximums at cxq ^ 
corresponding to interrupted regime when transitions between stable and metastable sliding friction 
are possible. 

The P a (a) dependencies in fig. [2] are plotted in log-log coordinates. It is seen that for curves 
DF and SS distribution takes the power-law form. Such regime corresponds to values ff< 1 and 



I a , I e <C It, at which (16) is written as follows 



P a {a) = a~ a V(a) 



(20) 



where V(a) is defined by: 



V(a) 



Z- l g- l I T 1/2 (l + a 2a ) x 



x exp 



r-i-; / l-ff[l-(2-T e )(l + (^ a ^-'i 

~ 1 rp C 



<T 9 



'1 + (a') 2a )- 2 (a') a 



da' 



(21) 



It is known that self-similar systems have a homogenous distribution [29]. Distribution (20) becomes 
homogenous at constant function (21 ). At small stress values the multiplier before exp is l+cr 2a — > 1. 



In fig.[3]the integration element is plotted without coefficient before integral. As can be seen from the 



figure integral in (21 ) has small value at a < 0.8, and when a exceeds certain value it rapidly begins 



to increase. According to structure of the equations (20), (21), integral gives the basic contribution 



in the resultant distribution (20), which becomes exponentially decreasing. Value a ~ 0.8 conforms 



with fig. [2] Thus, power-law distribution, typical for self-similar behavior, exists in limited range of 
stress values. Self-similar properties disappear when the stress exceeds the critical value. 



Stratonovich approach leads to first important difference from Ito form: in (20) distribution 
exponent is equal to —a, while in |T3j it is —2a. 




io 1 a 



Figure 3: Integrand of (21 ) at the same parameters as in fig. [I] and a = 0.5, 0.6, 0.7, 0.8, 0.9, 1.0 from 
bottom to top. 

3 Stress time series 



Euler method is used for numerical solution of the equation (12). Iterative procedure differs from 



analogous one used in |18| because Eq. (12) is the Stratonovich's SDE. To use common iterative 



procedure it is necessary to transform Stratonovich's SDE to the equivalent Ito's SDE. Taking into 



account the properties (13) for Eq. (12) one obtains the following form of the Ito's SDE 27,30 



da 







f a {<y) + VWi^VW) 



dt+^I a {a)dW{t). 



According to definition of the discrete analog of stochastic force differential dW(t) 



(22) 
^AtWi and 



(10), one can obtain iterative procedure for solution of the equation (22): 



a i+1 = a l + (fM) + a92< ^~ 1[ ^^P~ 2Ie] ) ^ + y/WfiiW*. (23) 

Solution of the equation runs over t e [0, T] time interval. At given numbers of iterations N 
(number of time series members) increment of time is defined as At = T/N. Force Wi has following 



peculiarities (cf. (13)): 



(Wi) = 0, {WiW v ) = 0, (Wf > -> 2. 



The Box-Muller model allows us to represent sufficiently stochastic force 31 

Wi = v / ^ v / -2hiri cos(2vrr 2 ), r n E (0, 1], 



(24) 



(25) 



9 



where, according to (24), /i 2 = 2 is the dispersion, Wi is the random number with properties (24). 



Pseudo-random numbers 7*1, r2 have uniform distribution and repeat themselves through periodical 



intervals. Effective potential (17) has minimums at positive and negative values of stress a. Thus 



while solving numerically Eq. (22) fluctuations cause transitions between states defined by the 



minimums. We can exclude negative part o < out of consideration since one-directional motion 
of upper moving surface is considered. This allows us to analyze further the behavior of |cr|(t). 
Typical realizations of the |cr|(t) for considered regimes are shown in fig. El Positive domains whith 




Figure 4: Stress time series |cr|(t), derived from equation (22) by numerical solution according to 
(23) at N = 10 4 ,t = 50, dt = 0.005. Regimes, that are shown in the plot, correspond to points in 



phase diagram (fig. [TJ . 

stress value close to zero is observed at dry friction regime (DF). There are random transitions 
between zero and non-zero stress values a at stick — slip regime (SS). Realizations of the SF 
and the MSF + SF regimes are visually similar. Therefore to detect the friction type one needs 
to apply additional analysis for probability density definition (see fig. |2|. Time series obtained in 
work [l~8] permit visual interpretation because corresponding phase diagrams were plotted in T e — I? 
coordinates, and a(t) were built at different values of friction surfaces temperature T e . At large 



10 



temperatures T e lubricant is totally melted, at small T e it is solid. Here, phase diagram is plotted 
at fixed value T e , therefore the time series related to different regimes are similar. They represent 
different friction regimes according to probability distributions shown in fig. [2j Note, plotting the 
phase diagram shown in fig. [TJ it is unreasonable to use large T e values because it is transformed 



into straight line (19) demarcating the SS and the SF friction regimes, and other domains are 



eliminated. For comparison, realizations of |cr|(t) are shown in fig. p^at the same parameters as in 
fig. |4| but at T e = 4. It can be defined visually that shown dependencies are in accordance with the 




Figure 5: Stress time series \cr\(t) corresponding to fig. plat T e = 4. 



SS and the SF regimes. There are transitions between zero and non-zero stresses values in the SS 
regime and in the SF regime always o > 0. However, the aim of this work is the analysis of the self- 
similar behavior, and we study all possible friction regimes. Therefore we use dependencies shown 
in fig. |4j In fig. [6] spectrum of the stress oscillations is shown, obtained by the fast Fourier transform 
algorithm (FFT) [31] at the SS regime time series analysis presented in fig. 111 Corresponding time 



series are obtained by iterative procedure (23) at iV = 2-10 ,t = 10 , dt = 0.005. It is evident 



that signal power in the spectrum is decreased with increase in frequency. White line described by 
relationship S p (f) oc l// 1 ' 8 is the spectrum approximation, i.e., power is inversely proportional to 
frequency. Thus there are different time correlations in the system, in contrast to the white noise 
for which S p (f) = const. For all considered regimes spectrums S p (f) have similar form, and for all 
cases S p (f) oc l// 1 ' 8 . Thus system considered on the basis of equation d9J) transforms white noise 
generators inherent in almost in all physical models into color noise with non-zero correlation time. 
Such behavior was observed experimentally (32). 



11 




Figure 6: Oscillation spectrum S p (f) corresponding to the SS regime parameters shown in fig. m 
White line is fixed by equation S p (f) oc l// 1,8 . Power S p is measured in conventional units. 

4 Multifractal detrended fluctuation analysis (MF-DFA) 

Multifractal analysis allows to calculate numerically the basic multifractal characteristics 1 33 1 de- 
scribing the self-similar systems. This method was proposed and developed by J. Kantelhardt et 



al. 34 , and it is widely used in many scientific fields for different time series analysis (meteorol 
o 



gy [35] , medicine 36 , economy 37 , and others [38}]40] ). Here we cite from original work [34] main 
statements of the method (for full description see [34]): Supposing that %k is a series of length N, 
to provide the multifractal detrended fluctuation analysis, one must follow next five steps: 



Step 1: Determine the "profile" 



Yd) 






<*>], 



1 N 



(26) 



Step 2: Divide the profile Y(i) into N s = int(Nfs) non-overlapping segments of equal lengths 
s. Since the length N of the series is often not a multiple of the considered time scale s, a 
short part at the end of the profile may remain. In order not to disregard this part of the 
series, the same procedure is repeated starting from the opposite end. Thereby, 2N S segments 
are obtained altogether. 

Step 3: Calculate the local trend for each of the 2N S segments by a least-square fit of the 
series. Then determine the fluctuation function 



F\v,s) = -j^{Y[(v-l)s + 



yS)Y 



(27) 



i=l 



12 



for each segments v, v = 1, ..., N s , and 

F>, s )= ] :J2{Y[N-(u-N s )s + l }- y^} 2 (28) 

s i=\ 

for v = N s + 1, ..., 2N S . Here, y v {i) is the fitting polynomial in segment v. Order of polynomial 
m selected with respect to the order of trend presenting in the series. Thus, polynomial of 
the m order can eliminate trend of order m — 1. 

• Step 4 '■ Average over all segments to obtain the q — th order fluctuation function: 

f -. 2N a ~\ V<2 

where, the index q can take any real value except zero. 

• Step 5: Determine the scaling behavior of the fluctuation function by analyzing the log-log 
plot of F q (s) versus s for each value of q. If the series Xi are long-range power-law correlated, 
F q {s) will depend on s, as a power law, 

F q ~ s h{q \ (30) 

where h(q) is the generalized Hurst exponent depending on q (note that h(q) at q = 2 is equal to 

classic Hurst exponent H |4l]). 

Function h(q) is connected with another classic multifractal scaling exponent r(g) |33~|[34] : 

r(q) = qh(q) - 1. (31) 

Self-similar behavior can be described by multifractal spectrum function f(a), connected with r(q) 
through Legendre transformation [33 j : 

a = r\q), f(a) =qa-r(q), (32) 

where a is the Holder exponent, and '"" denotes differentiation with respect to q. Using (31), 
can get directly related f(a) and h(q): 

a = h(q) + qh'(q), /(a) = q[a - %)] + 1. (33) 



we 



13 



Type of denoted dependencies characterizes time series behavior. Thus, constant value of h(q) = 
const and, correspondingly, linear increase in exponent r(q) denote monofractal series. Decrease in 
h(q) with q and nonlinear growth of r (q) are inherent in multifractal time series. Just one value of the 
Holder exponent a is characteristic for monofractal objects, and f{a) dependence presents a narrow 
peak. There is spectrum of f{ot) values in the case of multifractal series. However, in the case of 
monofractal series the numerical calculation will not give the only value of f(a), instead we will have 
a set of close values a, which, comparing with wider spectrums, approximately can be considered 
as monofractal issues. In general case, two types of multifractality can be distinguished for time 
series: (1) multifractality caused by broad probability density function of the series members, 
(2) multifractality caused by different time correlations between series members. Rearrangement of 
series components in a random order does not lead to the elimination multifractality in the first case. 
In the second case randomization leads to disappearance of the correlations. Since the reason for 
multifractality vanishes the series is transformed to monofractal. If both reasons of multifractality 
are inherent in series, the corresponding mixed series is characterized by weaker multifractality than 



initial one 34 . Thus, analyzing the mixed series using the method |34|, it is obviously possible to 
define the reason for multifractality and presence of time correlations. 

In general case, for time series two types of multifractality can be distinguished [34|: (1) multi- 
fractality caused by broad probability density function of the members of the series, (2) multifrac- 
tality caused by different range time correlations between series members. To define the reason for 
multifractality and presence of time correlations, one must apply shuffling procedure that consist 
in rearranging components of series in a casual order, and then compare corresponded spectrums 
f{a) of original and shuffled series. Thus, after shuffling of the series with multifractality of type 
(1) corresponded f(a) function (multifractality) must not changed, because the probability den- 
sity remains the same. In the second case hashing leads to disappearance of available correlations, 
and since the reason for multifractality thus vanishes, such series is transformed to monofractal. 
If for a series both reasons of multifractality are inherent in, the corresponding mixed series are 
characterized by weaker multifractality than initial series (34). 

Using this method we analyze stress time series \a\(t) shown in fig. Ill Figure [7] illustrates typical 
form of the F q (s) dependence, plotted in log-log coordinates, at some q values for time series related 
to the DF regime. From this figure we can see linear dependence on all set of s values, that is 
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Figure 7: Dependence of fluctuation function of F q (s) for the parameters of fig. Ill and the DF 



regime. 



typical for all series analyzed in the current work. 

This allows us to calculate precisely the Hurst exponent h(q) according to the scaling equation 



(30). We select domain 50 < s < 500 for calculation of the multifractal characteristics, where 



dependence of F q (s) has linear form. 

For time series shown in fig. Hat iV = 10 5 ,t = 10 3 , dt = 0.01 we can calculate h(q), r(q), and 
f(a). From fig. [8] it is seen that the strongest multifractality is exhibited by series related to the 
DF regime, then SS follows, and for series related to the MSF + SF and the SF regimes the weaker 
dependence h on q is characteristic that corresponds to monofractal behavior. Strong multifractality 
for the DF regime can be explained by power-law probability density function at small stresses, 
this is inherent in self-similar systems. In the SS regime multifractality is weaker, since probability 
density function of the Langevin equation solution also has non-zero maximum. There is possibility 
of system transition into state defined by this maximum, related to lubricant melting when it loses 
self-similar properties and settles into stationary regime of the sliding friction. But, existence of 
two maximums of -P a (c") corresponds to the stick — slip regime, and inverse transitions to solid-like 
structure may occur, and system returns to self-similarity 

For the MSF + SF and the SF regimes system does not exhibit multifractality, because of 
non-power-law probability density function. 

Peculiarity of the results shown in fig. [8] is that the multifractal characteristic for different 
friction regimes has close values in the q > domain and main difference observed in the area 



15 




10 



-5 



-10 



-15 



- 


i j# 


- 


«#^*^***2^ 


*dw 


a-DF 




o-SS 
*-SF 

o-MSF+SF 




Figure 8: Multifractal characteristics h(q), r(q), and /(a) corresponding to the parameters of time 



series in fig. U\ Set of curves 1 is related to series, derived directly from procedure (23), and set 2 
is related to analogous shuffled series. 



where q less than zero. In that range of q values MF-DFA takes into account small fluctuations in 
the time series, so it may seems that curves shown in fig. [8] may differs due to the peculiarity of the 
numerical realization of the MF-DFA procedure, and not due to the different statistical properties of 
the series, corresponded to variant friction regimes as we mentioned above. To explain this situation 
we must note that such results are typical for series with power-law distribution function. Thus, 
for uncorrelated multifractal series with distribution function 



P(x) = ax 



>+i) 



(34) 
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where a > and 1 < x < oo, corresponded Hurst exponent h(q) can be determined analytically |34|, 

and defines as: 

1/q (q > a) 

This dependencies shown in fig. [9| Latter curves have the same topology as results for stress 



h(q) 



(35) 




Figure 9: Hurst exponent h(q) for series with power-law distribution (34) obtained from eq.(35) 



time series shown on fig. [8j namely that main difference between multifractal Hurst exponent 
corresponded to series with different power-law distributions observed in q < domain. This is 
explicit result, characterful for time series with power-low distribution, and we suppose that results 
of our calculations of multifractal characteristics for stress time series are precision and correct. We 
also note that for stress time series magnitude of the corresponded h(q) function not determined by 



power exponent like in eq.(35), because in that case we have power-law distribution only for small 



stress values (for DF regime approximately a < 1), and for larger a it breaks (see fig. pi). 

According to above mentioned, we can conclude that in this case multifractality is caused by 
power-law distribution function. To detect different time correlations which may be present in the 
system, we need to shuffle series and then, again calculate multifractal characteristics. In figure 10 
spectrum of the stress oscillations related to shuffled series, described by fig. [6j is shown. White line 
is the approximation of spectrum and can be described by relationship S p (f) oc l// 0017 , i.e., power 
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f,Hz 



Figure 10: Oscillation spectrum S p (f), related to shuffled series described by fig. p| White line 
corresponds to relationship S p (f) oc 1//° . Power S p is measured in conventional units. 

almost does not depend on frequency. It means that correlations in system are disappear. Given 
spectrum is related to white noise, S p (f) = const. Thus, shuffling of the series causes elimination 
of correlations. But, since while shuffling the time series, neither addition nor substraction of the 
series members are performed, distribution function stays the same. 

Set of curves 2 in fig. [8] corresponds to shuffled series while set 1 is related to original ones. As 
we see h(q) is a straight line h = 0.5, the spectrum function f(a) presents a narrow peak with 
small width, and r(q) is a straight line with constant slope. Denoted peculiarities are related to 
monofractal system, and value h = 0.5 corresponds to uncorrelated series. Thus, for considered 
system multifractality is caused by power-law distribution function and by different correlations. If 
power- law dependence of P a {(x) is broken, or correlations vanish, the multifractality is eliminated. 



5 Conclusion 



Using the homogenous rheological model ultrathin lubricant film melting has been investigated. 
Basic parameters are shear stress and strain, as well as temperature of lubricant. Four regimes of 
lubricant behavior, characterized by different sets of maximums of stresses distribution function, 
have been found. Stress time series have been obtained for each regime by numerical modeling of the 
Langevin equation, and it has been shown that at specific parameters time series are multifractal. 
All basic multifractal characteristics have been calculated, and it has been shown that multifractality 
is caused by different time correlation and also by power-law distribution in the limited range of the 
stress values. When temperature noise intensity is much larger than intensities of stress and strain 
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noises a power-law distribution can be observed. According to above examination, multifractal time 
series are realized only for the dry friction (DF) and stick — slip (SS) domains, since only for this 
regimes power-law distribution is observed. 
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